Reconstruction of Binary Functions and Shapes 
from Incomplete Frequency Information 

Yu Mao 

Institute for Mathematics and Its Applications, University of Minnesota. 



(N 

O 
(N 



(N 



C/2 



> 

o 
o 



Abstract — The characterization of a binary function by partial 
frequency information is considered. We show that it is possible to 
reconstruct binary signals from incomplete frequency measure- 
ments via the solution of a simple linear optimization problem. 
We further prove that if a binary function is spatially structured 
(e.g. a general black-white image or an indicator function of a 
shape), then it can be recovered from very few low frequency 
measurements in general. These results would lead to efficient 
methods of sensing, characterizing and recovering a binary signal 
or a shape as well as other applications like deconvolution of 
binary functions blurred by a low-pass filter. Numerical results 
are provided to demonstrate the theoretical arguments. 



I. Introduction 

This paper discusses the reconstruction of the binary signals. 
Binary signals appear in a variety of applications like shape 
processing, bar code and handwriting recognition, obstacle 
detection, image segmentation; see e.g. ll24l . ll28ll . 1201 . ifTSll . 
[[Tl, 1 29 1 and many others. 

One of the major difficulties in the reconstruction of bi- 
nary functions is that the binary constraint is non-convex. 
Optimization with a binary constraint is often approached by 
means of the double- well potential or other nonlinear schemes. 
In this paper we demonstrate that binary functions can be 
reconstructed exactly via a simple convex optimization when 
only partial frequency information is available (e.g. when the 
signal is blurred by a low-pass filter). 

A. Main results 

Let uo be a binary function, i.e. uo{x) G {0,1}, Vx. 
Let T be the Fourier transform and S the selecting operator 
corresponding to the incomplete measurements b. Our goal is 
to recover uq from b = STuq, which is an underdetermined 
problem. The main contribution of this work is showing that 
under certain conditions, uq can be exactly reconstructed by 
solving the convex relaxed optimization problem 

find u s.t. STu = 6, < u < 1. 

At first glance, this may seem a little bit surprising, as it is 

not even obvious that the solution of the problem is unique. 

However, in this work we prove that in many cases, the 

solution is unique and is equal to uq. 

• When a binary signal is spatially structured, i.e. the 
Is and Os are clustered (e.g. in a binary image or 
as an indicator function of a shape), with very few 
low frequency measurements taken, the solution of this 
convex optimization problem is deterministically unique 



and equals to the original binary signal. For a detailed 
statement, see theorem HOI 

• If a binary signal has no spatial structure, for example 
if the Is and Os appear randomly, we show that this 
relaxation works with overwhelming probability when the 
number of the measurements is more than a half of the 
size of the signal, and the probability tends to 1 as the size 
of the signal increases to infinity. For a detailed statement, 
see theorem [II. 141 

• We also propose a very efficient algorithm designed 
for this convex problem (see algorithm [T}. Numerical 
experiments are presented in section |IV] 

B. Related works 

The idea that under certain circumstances, the binary con- 
straint can be automatically satisfied by imposing a convex 
relaxation, in particular the box constraint < < 1, is 
not new. For example, when solving the image segmentation, 
multi-label and many other problems based on the total 
variation model (see e.g. 0, O, l|29l ), people have noticed 
that although the original problem is non-convex, the global 
minimizer can be obtained by solving the relaxed convex 
problem and the solution will be automatically (almost) binary. 
However, this approach works only because of the special 
structure of the variational model. The theoretical analysis 
strongly depends on the coarea formula for total variation. 

In the contexts of regression and approximation, it is has 
been known for a long time that in an Loo regression (some- 
times referred to as Chebyshev or minimax regression where 
the penalty function is given by the Loo norm) or a deadzone- 
linear penalty regression (where the penalty function is given 
by the deadzone function (| • | — a)+), the distribution of the 
residual of the regression will concentrate at the boundary 
of the feasible domain (see e.g. 0|, Chapter 6). When the 
feasible domain is interval [0, 1], a function with many values 
right at the boundary is nothing but a binary signal. In fact, 
in section III-AI we will show that our convex relaxation of 
the problem can be equivalently reformulated as an Loo or 
deadzone penalty minimization problem. 

The idea of recovering a signal from the partial frequency 
measurements is often used for compressed sensing [5 |, [6|, 
1 12 1, which takes advantage of the prior assumption on the 
sparsity of the signal and reconstructs the signal via Li 
minimization. However, the current work is substantially dif- 
ferent from compressed sensing. Although structured binary 
functions are a special case of piecewise constant functions 
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whose derivative is sparse, the condition of being binary 
is actually stronger than simply being piecewise constant, 
therefore stronger results can be expected. Indeed, none of 
the results given in this work can be deduced from the 
standard compressed sensing theory directly, and some of them 
are of a very different nature. For example, in compressed 
sensing the frequency measurements should be taken randomly 
to guarantee the restricted isometry property [5 |, while in 
the reconstruction of the structured binary function, the low 
frequency measurements actually play a more important role 
than the high frequency measurements as discussed below. 
On the other hand, many major results given in this paper 
are deterministic, while results in the compressed sensing 
literature are often intrinsically stochastic. 

The present research is also related to a seminal work on the 
reconstruction of signals from partial frequency information 
03 1, where the spatial structures and patterns in both time and 
frequency domain are used to guarantee the uniqueness of the 
signal reconstruction. In the very recent research [im . [lT4ll . the 
authors showed that a random binary signal can be recovered 
with certain probability by means of the relaxed box constraint. 
In their work the major mathematical tool is the delicate 
geometric face-counting of random poly topes. We also get a 
basically similar result in section III-Fl but from a different 
approach. In ll36ll . the authors defined the degrees of freedom 
contained in a sparse or piecewise polynomial signal as the 
rate of innovation of the signal. Then they showed that the 
quantity of the samplings needed to recover the signal equals 
to the rate of innovation. However, the mathematics behind 
their theory is substantially different than ours. Moreover, 
their method requires a certain pattern of sampling and the 
reconstruction involves a factorization of polynomial. In [41 
the authors showed that if an underdetermined system admits a 
very sparse nonnegative solution and the matrix has a row-span 
intersecting the positive orthant, the solution is actually unique. 
In 06] the author proved that a sparse nonnegative can be 
reconstructed as the unique solution of a linear programming 
problem, where the corresponding matrix is the submatrix of 
a Fourier matrix consisting of its top rows. We will further 
discuss the relationship between these two works and ours in 
section III-CI 

C. Notations and conventions 

In this paper all signals are assumed to have periodic 
boundary condition. The /i-dimensional discrete signals are 
defined on {1,...,A^}^ where N is assumed to be even. 
Here for simplicity we assume that the domain is equilong 
along each dimension. The /i-dimensional continuous signals 
are defined on = [0, 1]^ where the two endpoints and 1 
are identified due to the periodic boundary condition. 

We use to denote the Fourier transform, both in the 
discrete and continuous periodic cases. 

When we talk about the discrete Fourier transform, we use 
the following convention: 



where {a/e} are the Fourier coefficients defined on a symmetric 
support (N/2 is treated as same as —N/2). A smaller 
corresponds to a lower frequency. Since u is always real, {ak} 
satisfies a-k = Notice that = N^T~^. 

We use S to denote the selecting operator. 5 is a diagonal 
matrix where the selected positions have value 1 and others 
are 0. 

D. Contents 

The paper is organized as follows. Section HI] discusses 
the theoretical results. In section [IIIl an algorithm to solve 
the convex problem is proposed. Numerical experiments are 
shown in section [Vl and conclusion is given in section |Vl 
To make the main text more concise, we put all proofs into 
Appendix except those theorems and corollaries immediately 
deduced from the discussion in the context. 

II. Theory 

A. General reconstruction theory 

Suppose is a discrete binary signal, i.e. uq{x) G 
{0, 1}, Vx. Consider a linear system Auq = h where A = ST, 
T is the Fourier transform and S is the selecting operator. 
The meaning of this system is clear: some partial frequency 
information of the binary signal is given, and we want to 
reconstruct iio from the incomplete measurements. This leads 
to the following problem (Pq): 

Po : find u s.t. Au = b, u{x) e {0, 1}. (1) 

The problem (Pq) is non-convex due to the binary condition, 
and the following convex problem is the tight relaxation of 

Pi : find u s.t. Au = b, 0<u<l. (2) 

We want to show that (Pi) can be used to recover uq 
under certain conditions. The following theorem specifies the 
conditions guaranteeing that this relaxation is exact. 

Theorem II.l. Assume uq is a binary solution of Auq = b. 
There exists no nonzero v e {Av = 0} such that 

v{x) < 0, when uo{x) — 1 
v{x) > 0, when uo{x) = 

if and only if uq is the unique solution of (Pi), i.e. solving 
(Pi) recovers uq. 

If the size of the signal is N, then the criteria © on v 
determines an orthant in depending on uq. We denote this 
orthant as O^o • The theorem tells us that as long as the kernel 
space of A intersects O^o at nowhere but the origin, solving 
(Pi) is enough to recover uq. 

This condition is 'negative', i.e. it requires the nonexistence 
of such a vector v. The following statement, sometimes 
referred to as the Gordan-Stiemke theorem of the alternative, 
will lead to a 'positive' criteria. 
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Lemma II.2 (Alternative Theorem, see e.g. [T|). One and 

only one of the two following problems is feasible: (1). Find 
0^v>0 s.t. Av = 0; (2). Find v = A^r] s.t. v >0. 

Geometrically, this theorem says that if P is a subspace 
in R^, O is the first orthant, then either Pp|0 = {0} or 
Pi int(O) = 0, but not both. The statement considers the 
first orthant only, but obviously it is true for any other given 
orthant. Apply this lemma on theorem III. II we immediately 
get the 'positive' version of the criteria: 

Theorem II.3. Assume uq is a binary solution of Auq = 
b.There exists v = A^r] such that 

v{x) < 0, when uo{x) = 1 
v{x) > 0, when uo{x) = 

if and only if uq is the unique solution of (Pi), i.e. solving 
(Pi) recovers uq. 

Unfortunately, there is no explicit formula to determine if 
an arbitrary subspace passes through a given orthant. Indeed, it 
is equivalent with any general linear programming feasibility 
problem and thus has no closed-form solution. However, for 
some special cases we can still give deterministic or stochastic 
results, as we will explain in the following subsections. 

We want to remark that there are many alternative linear 
programming problems that can recover the signal as well. In 
fact, assume J{u) is a convex function on u satisfying 

J{u) < J{v),yu e [0,1]^,^^ [0,1]^. (5) 

It is easy to see that if uq is a binary solution of Auq = b, 
then Uq is a unique solution of (Pi ) implies that uq is a unique 
solution of the following convex problem: 

min J{u) s.t. Au = b. (6) 

u 

Therefore solving ^ can also recover uq under the condition 
in theorem UlTI or III. 3l There are many simple functions satis- 
fying (O. One of the simplest examples is J{u) = \\2u — l\\oo. 
Another example is the deadzone penalty J{u) = J{{\2u — 
1| — 1)+) where J is any convex function with J(0) = and 
J{v) > for V 7^ 0, e.g. J{v) = \\v\\p for p > 1. 

B. Reconstruction of the ID binary signals from the low 
frequency measurements 

So far the discussion has used only the fact that the signal 
to be reconstructed is binary. In most practical applications, 
the signal is often not only binary, but also structured, i.e. the 
Is and Os are spatially clustered. This property could help us 
reconstruct the signal. 

Let us consider the ID case first. Assume uo{x) is a 
periodic discrete binary signal defined on {1, . . . , A/"}. Since 
we are considering the structured signal, uo{x) consists of 
many intervals with constant value 1 or 0. If the first and last 
intervals are with the same value, we treat them as one merged 
interval under the periodic boundary condition. Therefore, the 
total number of intervals is always even, thus uo{x) can be 



represented as 

2d 

{Ij} is a partition of {1, . . . , N} where each Ij is a consec- 
utive interval. 

From theorem III. 31 uq can be recovered from (Pi ) if and 
only if there exists v = A^r] such that 

v{x) < in if = 1 
v{x) > in Ij if = 

We want to show that this condition is always satisfied 
for certain types of A. Recall that when partial frequency 
information is given, A = SF where 5* is a sampling operator 
that corresponds to the known frequencies. If v = A^r] = 
[St]), then v is a band-limit signal whose spectrum can be 
represented as Sr], i.e. it is located inside the known frequen- 
cies. Therefore, the above condition means that the relaxation 
method is valid as long as we can use only those known 
frequencies to construct a band-limit signal that satisfies ©. 
Since ^ describes the zero-crossing position of it imposes 
a constraint on the spectrum of v, and therefore on S. 

The relationship between the zero-crossings of a signal 
and its spectrum information is not a new problem in signal 
processing; readers are referred to [25 , [30 , [21], [35 for 
some classic theories. The following result is natural from the 
perspective of trigonometric interpolation: 

Lemma IL4. Let T = [0, 1] where and 1 are identified, i.e. 
T = S'^ = {z : \z\ = 1}. Given 2n points on T who define 2n 
intervals on T, there exists a real trigonometric polynomial, 
whose spectrum is limited in [— n, n], vanishing only at those 
points and changes signs alternatively on those intervals. 

This conclusion, combined with theorem III. 31 leads to the 
following deterministic result which states that the number of 
low frequency measurements we need to reconstruct the binary 
function is basically the number of the jumps contained in the 
signal, no matter how large the signal is. 

Theorem IL5. If uq{x) is a 1-D binary signal that can be 
represented as in d?]) with 2d consecutive intervals of ones and 
zeros, then by knowing the Fourier coefficients {ak} for \k\ < 
d, we can recover uq through the convex problem (Pi). (Notice 
that uo{x) G ]R,Vx implies ak = cT-k^Mk, so essentially we 
only need to know {ak} for < k < d.) This result is optimal, 
i.e. precise reconstruction via solving (Pi) is impossible if 
knowing even less. 

Although theorem III. 51 only holds when the lowest fre- 
quency information is given, it is still very useful, because 
in many practical problems the low frequency measurements 
are far easier to obtain than the high frequency measurements. 
The deconvolution problem with a low-pass filter kernel, for 
example, can be treated as reconstruction from the lowest 
frequency information. 

Heuristically, theorem llL5l can be understood as follows: if 
the low frequency measurements are given, then the permitted 
perturbation can be with higher frequencies only and thus 
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Strongly oscillating around zero. Therefore, by controlling the 
lower and upper bounds of the signal as in (Pi), the oscillating 
perturbation would be eliminated, and thus the solution is 
uniquely determined. 

C. Discussions and generalizations 

First, it is easy to see that theorem III. 51 can be directly 
extended to the cosine transform as well, due to the fact that 
the cosine transform of a signal is nothing but the Fourier 
transform of the even extension of the signal. 

Corollary II.6. If uq{x) is a 1-D binary signal that can be 
represented as in ^ with 2d consecutive intervals of ones 
and zeros, then by knowing the discrete cosine transform 
coefficients {au} for < k < 2d, we can recover uq through 
the convex relaxation (Pi). 

We also mention that if the signal is only bounded from 
one side, i.e. instead of knowing that the signal is binary, we 
know the signal is nonnegative, then a similar argument would 
lead to a theorem concerning the reconstruction of the sparse 
nonnegative signals. Indeed, using theorem III. II and lemma 
III. 21 we can obtain the following theorem (the proof is similar 
hence omitted): 

Theorem II.7. If uq > is supported on K = {x : uo{x) ^ 
0}, then uo is the unique solution of Au = b, u > if and 
only if there exists v = A^r] such that v\k = 0, v\k^ > 0. 

Let A = ST and S also select the low frequency mea- 
surements, then theorem III. 71 implies the following theorem 
(thanks to lemma llL4l as well): 

Theorem II.8. If uq{x) is a 1-D nonnegative sparse signal 
supported on K = {x : uo{x) = 0} with \K\ = d, then by 
knowing the Fourier coefficients {ak} for \k\ < d, we can 
recover uq through the convex problem Au = b, u > 0. 

This result is closely related with the theorems proved in 
lUl which said that a nonnegative solution of a linear system 
is unique if the solution is sparse enough and the matrix 
has a row-span intersecting the positive orthant. It is worth 
mentioning that the quantity of the needed low frequency 
measurements in this case approximately equals two times the 
quantity of the 'spikes' of uq (other than the number of jumps 
of Uo in theorem HO) . This coincides the observation in [36 | 
that the degree of freedom of a sparse signal is 2d (for each 
spike there is one degree for position and one for amplitude), 
and thus 2d -\- 1 measurements are in principal enough. A 
similar observation has been given in [[T6l as well. 

Our result may further be generalized to bases other than 
the trigonometric functions. Indeed, the duality of lemma llL4l 
tells us that a signal without lower frequency components must 
have many sign changes, which is some times referred to as the 
Sturm-Hurwitz theorem [35 1. This observation plays a critical 
role here. This property can be extended to other basis that has 
similar oscillating pattern, such as some wavelet bases or the 
eigenfunctions of the regular Sturm-Liouville problems [17l. 
However, generalization along this line is beyond the scope of 
this paper. 



D. Reconstruction of the 2D binary signals from the low 
frequency measurements 

The multidimensional case is more complicated than the 
ID case due to the following several reasons. There is no 
fundamental algebraic theorem for multivariable polynomials. 
Moreover, the Sturm-Hurwitz theorem that describes the zero- 
crossings of function with a spectrum gap does not exist 
in higher dimensions. Finally, in ID the complexity of a 
binary function can be simply characterized by the number 
of jumps as in theorem III. 51 while in higher dimensions, a 
binary function may have only one connected component but 
still have a very complicated jump set. 

Since theorem III. 31 is still valid in the multidimensional 
case, a multidimensional binary signal can be reconstructed 
by the lower frequency measurements as long as the jump set 
of the binary signal is the zero levelset of a low frequency 
function. Unfortunately, to the author's knowledge, no criteria 
has been known to determine if a given shape can be realized 
as the zero levelset of a function with only lower frequency 
components. In (13, d, O, (33, O, O, E3, some 
results concerning the relationship between the levelset of a 
function and its Fourier transform are shown. In [26] it has 
been proved that using the continuous Fourier transform, a 
function with a given levelset curve can be approximated 
to any degree of accuracy by a band-limited function with 
given spectrum support. However, this result is barely useful in 
practice because it requires virtually infinitely high resolution 
in the frequency domain. 

Heuristically, if a function has only lower frequency com- 
ponents, we can imagine that its levelset would not be too 
complicated. Here we give a way to measure this complexity. 
The basic idea is that since in ID case the complexity of a 
binary signal is determined by the number of jumps inside the 
signal, in 2D case we can define an 'average number of zero- 
crossings' as illustrated in Fig. [T] The following discussion 
can be naturally extended to higher dimensional cases. 

Define = [0, 1] x [0, 1] with opposite boundaries identi- 
fied, that is, T2 ^ R7z2. For e (-7r/4, 7r/4], define 

Ls,e{t) = (t, s + nan i9) mod 1, s e [0,1], t G [0,1] 

to be a grating along angle 6 starting from the left edge of the 
square (see Fig. [T] left). For G (7r/4, 37r/4], similarly define 

Ls,e{t) = (s + t cot i9,t) mod 1, s e [0,1], t G [0,1] 




Fig. 1. Two gratings on with a binary function. 
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to be the grating along angle starting from the bottom edge 
of the square (see Fig. [T] right). Assume is a binary function 
whose jump set consists of analytic curves. For given 9 and 
s the line segment Lg^e will intersect the jump set of u finite 
times. We denote this number as #Ls 51 and define the average 
of #Ls^e over s as 

Ke=cosO / ij^Ls^eds 
Jo 

for e (-7r/4,7r/4] and 

Ke = sm6 / #Ls^eds 
Jo 

for e (7r/4, 37r/4]. The presence of the multiplier cosO 
and sinO is due to the fact that for gratings with different 
angle 0, s is not an equilong variable. The distance between 
Ls^e and the origin is more intrinsic which equals scosO for 
e (-7r/4,7r/4] and ssinO for G (7r/4, 37r/4]. Kq is 
called the average directional number of zero -crossings in this 
paper which essentially describes the average quantity of sign 
changes along the direction 6. It is easy to see the connection 
of this quantity and the number of jumps for a ID binary 
signal. Indeed, by the Cauchy-Crofton formula, J KedO is 
nothing but two times the perimeter of the shape, i.e. the total 
variation of uq, while in the ID case the number of jumps also 
equals to the total variation of the binary signal. In the next 
theorem we will show that Kg in some sense characterizes the 
complexity of the shape. 

Theorem II.9. Assume u{x^y) is a 2D binary function with 
analytic jump curve and the average directional number 
of zero -crossings of u along the angle 9 is denoted as 
Kg. If there exists a band-limited real function v{x^y) = 
Eo-,/c)G^^^ife^^'''^^''^^^^ defined on T^, where = {{j,k) : 
p + k'^ < d}, such that the jump set ofu{x^ y) corresponds 
to the zero levelset of v{x^y), then Kq < 2d,\/0. 

Corollary 11.10. Assume uo{x^y) is a discrete 2D binary 
function defined on {1,...,A/'}^ and u{x^y) is a binary 
function defined on the continuous domain with ana- 
lytic jump curves such that u{x/N^y/N) = uo{x^y) for 
(x^y) e {1,...,A^}^, and denote the average directional 
number of zero-crossings ofu along the angle by Kq. If the 
reconstruction ofuo by linear programming problem {Pi) from 
low frequency measurements in ft = {(j, k) : y/p^rk^ < d} 
is exact, then d> ^ max^i Kq. 

The meaning of theorem III. 91 and corollary III. 101 is clear: 
the average directional number of zero-crossings of the levelset 
of a band-limited function is bounded by the diameter of the 
support of the spectrum. If we denote the jump set of u{x^ y) 
as r and the perimeter |r|, then by Cauchy-Crofton formula, 
the condition d>\ max6i Kq further implies d> :^ |r|, which 
can be seen as a natural generalization of the ID case (see 
lemma III.4I and theorem III. 51) . However, unlike the ID case, 
this theorem just gives the necessary condition, not a sufficient 
one. 



E. Reconstruction of binary signal from arbitrary frequency 
measurements 

If S is an arbitrary frequency selector, not necessarily se- 
lecting the lowest frequencies, it is not easy to give a sufficient 
and necessary condition to determine if the reconstruction is 
possible since there is no way to quantify the zero-crossings 
simply from the irregular support of the spectrum. In [! 23l 
the authors show that given the support of the spectrum of 
a trigonometric polynomial, the size of the largest non-zero 
circular region of the polynomial is bounded. They proved the 
following theorem: 

Theorem 11.11. Let ^ S C Z"^ be a finite set s.t. S = -S. 

Let v{x) = "^Zkes c^e^^*^^'^^ be a real valued trigonometric 
polynomial on T^, then v{x) has at least one zero in any 
closed ball of diameter 

This theorem indicated that if the known frequencies have 
an arbitrary support, then the binary functions can be recovered 
from (Pi) if it contains a constant block large enough. 
However, the bound given in theorem III. Ill is rather loose. 

It is worth mentioning that the conclusion of theorem III. Ill 
tells us the 'importance' of each frequency band is roughly 
determined by the reciprocal of the frequency. That is to say, 
knowing lower frequency measurements is more important for 
reconstruction of the binary signals than knowing the high 
frequency measurements. This coincides with the intuition 
we learn from theorem III. 5 1 and differs from the case of 
sparse reconstruction as in the compressed sensing problems, 
where the measurements should be spread out in the frequency 
domain as much as possible. 

F. Reconstruction of random binary signal 

If a binary function is random, i.e. the orthant O^o is 
randomly chosen, there is no deterministic way to guarantee 
if a certain subspace passes through it, but the probability can 
be estimated. From now on we denote the kernel of A, the 
image of A'^ by I a and the rank of A by Ka, I a and r 
respectively, then dim (i^T^) = N — r, dim(/^) = r. We say 
an r-dimensional linear subspace is in general position if the 
projections of any r axes of onto the subspace are linearly 
independent, and we say A is in general position if Ka is 
in general position. The following result has been known by 
mathematicians at least as far back as the 1950s (see |8| for 
a brief review). It says that any r-dimensional subspace in 
in general position will pass through a fixed number of 
orthants of R^: 

Lemma 11.12 (see e.g. fSl). Any r-dimensional subspace in 
R^ in general position passes through 2 ^^Zq {^7^) orthants 
ofR^. 

We denote Pr^N = Yll=o (T)/^^' which is nothing but 
the cumulative distribution of the function of the standard 
binomial distribution with p = \- For a random binary signal 
uq, we say uq has no 'preference' on orthants if for any two 
orthants Oi and O2, 

Prob(Ono = = Prob(Ono = ±02), (9) 
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Fig. 2. Pr,N as a function of r/N. 



then since there are 2^ orthants in total, we have 



2iv 



^r-l,Ar-l 

(10) 



According to theorem HOI this is equivalent to saying: 

Theorem 11.13. If uq is a random binary signal with size 
N without preference on orthants, then given a matrix A in 
general position with rank r, the probability that uq can be 
recovered from linear problem (Pi) is Pr-i^N-i- 

It is well known that Pr^N can be approximated by 



1 r 



is the cumulative distribution function of the normal distribu- 
tion. By Hoeffding's inequality, the tail of Pr^N is bounded 
by 



Pr .N ^ 1 ^ 



when r < N/2 



. (11) 



f_{2r-Nf\ 

"r.N r:: J- - f exp i- ^'^''2N^ ) whcn r > N/2 
Therefore, if r/N ^ pas N ^ 00, then 

Pr,N <lexp(^-{p- N^ ^0 when p < 1/2 

PrAT > 1 - ^exp (- (p- ^)^7V) ^ 1 when p > 1/2 ' 

^ ^ (12) 

which is illustrated in Fig. [21 

The discussion above can be summarized by the following 
theorem, which basically says that if the number of measure- 
ments are more than a half of the size of the signal, the 
probability that the convex relaxation is exact will tend to 
1 as the size of the signal goes to infinity. 

Theorem 11.14. If uq is a random binary signal with size N 
without preference on orthants, A is a matrix in general posi- 
tion with rank r, when N is large, the probability that uq can 
be recovered from linear problem (Pi) can be approximated by 
^ (^ ^^N-i^ ^ w/z^r^ ^ is the cumulative distribution function 
of the normal distribution. If p > 1/2 as N ^ 00, 



then Uq can be recovered from (Pi) with overwhelming 



probability at least 1 



Ig c(Ar 1) y^li^j.^ g . 



III. Solving the Optimization Problem 

As stated in section UTAI there are many convex models that 
can recover the binary signals. Since the measurements might 
be noisy in practice, we choose to reconstruct the signals via 
the following optimization problem: 



min \\Au-bf s.t. 0<u<l. 



(13) 



First we discuss the robustness of this model. Let is the 
true binary function and b = Auq is the clean measurement. 
Assume b is contaminated by noise e. For corrupted measure- 
ment b = b-\- e, we want to investigate if this model will still 
lead to the correct answer. Let 



B{u) 



1, Ui> 1/2 
0, Ui < 1/2 



be the thresholding operator that maps any function to its 
closest binary function. The following theorem shows that 
the model is robust to small perturbation. In section |IV] the 
numerical results will show that the more measurements are 
given, the more robust the reconstruction would be, which is 
not surprising. 

Theorem III.l. If uq is the unique solution of (Pi), b = 
b-\- e is the corrupted measurement, u is the minimizer of the 
optimization problem 

min \\Au-b\\l s.t. 0<ii<l, 



where h > is a small amount 
Duo (^^^ details in the proof). 



then when \\e\\ < h{A,Ouo) 
depending only on A and 

B{u) Uq. 

Since ([T3]) is a standard bounded least square problem, it can 
be solved by many existing optimization algorithms. However, 
we propose an algorithm that is specifically developed for this 
problem. It will only utilize the discrete Fourier transform 
without explicitly storing and multiplying the matrix A, which 
is can be very large in practical problems and can make most 
out-of-the-box optimization packages very inefficient. 

Our algorithm will be based on the split Bregman method 
introduced in iflSl and modified in |34| for solving the 
non-negative least square problem. We replace (fT3l l with an 
equivalent problem 



min IIAii — s.t. 



p{d). 



where P{d) is defined component- wisely by 

P{d) = < 



1 d>l 
d < d <1 . 
d<0 



This constrained problem can be solved iteratively by 



mmd,u^\\Au-b^\\l^\\u-P{d) 



.k\\2 



- V 
b^ 



-P{d^^^) - 

b - Au^^^ 



The last two lines are called Bregman steps and can be 
understood as the gradient ascent steps in the augmented 
Lagrangian method. Theory on the convergence of this method 
can be found in ll22ll . ITTll , ifTSll . The first line can be solved 
exactly respectively on d and u, giving rise to the following 
iterations: 

^fe+l = 5^ ^ ^ _ 

(14) 

Here the first, third and last lines contain only trivial computa- 
tions. For the second line, we can notice that when A = SF, 
(AA^A + /)-i = J-i(ArA5' + /)-iJ' where N is the size of 
the signal. Since N\S + / is nothing but a diagonal matrix, 
the whole operator can be calculated efficiently and precisely. 
Therefore we get the following algorithm [T] Numerical results 
in the next section would show that this algorithm works very 
well. 

If the given information is not the partial Fourier mea- 
surements of the signal but a filtered signal, i.e. h = Au = 
KTu where i^^ is a filter in frequency domain, then this 
algorithm still works after a small modification. The only 
point that needs to be changed is that now (XAJ A + /)~^ = 

In some problems, the norm \Au — 6|p can also be precon- 
ditioned, e.g. to prevent the effect of the noise in high frequen- 
cies. We can minimize \^Au — = (^'^ ~ hY M{^Au — b) 
where M is a certain preconditioner in the frequency domain. 
Again, the algorithm still works without many modifications 
except the inverse operator becoming {XA^ MA + /)~^ now. 

Algorithm 1 The Split Bregman Algorithm for Solving (fT3l l 
Initialize: Let bo = b. Start from initial guess u = T~^b. 
while \\Au — ^o|P not small enough do 

d ^ P(u - v) 

u ^ {\A~^A + I)-^{\A~^b + P{d) + v) 

V ^ V ^ P{d) - u 

b ^ b^bo - Au 
end while 
u ^ B{u). 



IV. Numerical Results 

First we numerically show that binary signals can in general 
be characterized by very few frequency measurements. Fig. |3] 
shows several ID and 2D signals: 

• The first one is a ID binary signal. It contains 15 
constant intervals with value 1 and 15 intervals with 0, 
so by theorem III. 51 it can be fully determined by Fourier 
coefficients ak for \k\ < 16, no matter how long the 
signal actually is (the length of the signal shown here is 
400). 

• The second one is a binary image corresponding to a 
geometrical shape, the size is 200 x 200. The experiment 
shows that it is fully determined by Fourier coefficients 
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Fig. 3. Several binary signals: a random ID binary signal, a geometrical 
shape, a barcode and a handwriting image. 




Fig. 4. The signals after low-pass filtering and the reconstructions. 

a/c for \k\ < 5. The needed Fourier coefficients for 
characterization of the shape account for 0.2% of the total 
Fourier coefficients. 

• The third one is a barcode image, which can also be 
treated essentially as a ID signal. It contains 15 black 
bars and 15 white bars, so the Fourier coefficients needed 
are ak with < 16, where hi is the component of k 
along the horizontal dimension, no matter how large the 
image actually is (the width of the barcode shown here 
is 400). 

• The last one is an image with handwritten letters, with 
size 100 X 100. The experiment shows that it is fully 
determined by Fourier coefficients ak with \k\ < 10. The 
needed Fourier coefficients for the characterization of the 
image account for 3.17% of the total Fourier coefficients. 

To show that these binary signal can be recovered from 
the low frequency measurements, we filter them with a low- 
pass Gaussian kernel whose band corresponds to the partial 
frequencies. As long as the needed low frequency information 
is precisely given, an exact reconstruction would be available. 
However, the numerical experiments also show that the more 
measurements are known, the faster the reconstruction is, 
which means that it is easier for the algorithm to find the 
correct binary signal. Fig. |4] demonstrates the filtered signal 
and the reconstruction. The curves in Fig. |5] show that the 
reconstruction time decreases when more measurements are 
given. The x-axis measures the radius of the support of the 
given frequency information, while the axis measures the 
logarithm of computational time. 

To demonstrate that the reconstruction is robust, we now 
recover the signal from the low-pass filtered measurements 
with noise. It is well known that deblurring with large amounts 
of noise present is a difficult task. Our results show that 
even with very large amounts of noise and strong blurring, 
the results are still sensible. In Fig. |6l tests on a ID binary 
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Fig. 5. Time elapsed (in seconds) for the reconstruction versus the extra 
radius of the support of the given frequency information. 



noisy signal 




reconstructed signal 



signal with noisy measurements are demonstrated. When the 
input signal is very noisy, the positions of the bars in the 
reconstruction are not precisely equal but close to the original 
signal. The minor difference between the reconstruction and 
the true signal are highlighted by circles. Fig. |7] demonstrates 
how the number of miss-identified Os and Is changes with 
the number of measurements and noise level. We do the 
experiments for different levels of noise and different numbers 
of measurements respectively. For any given pair of fixed noise 
level and number of measurements, 1000 random tests are 
taken to get the average number of miss-identified Os and Is. 
In Fig. [SlfTOl each group shows a signal with a different level 
of noise. The parameters of the blurring kernel and the noise 
levels are given in the captions. The computational costs are 
also recorded. All computations are done in Matlab on a 2.8 
GHz Intel CPU. 



noisy signal 




reconstructed signal 



i9 rn n Qi 



noisy signal 



V. Conclusion 

The reconstruction of binary functions is difficult because 
of the nonconvex nature of the problem. In this work we 
proved that even with very few frequency measurements, 
binary functions can actually be reconstructed by solving a 
very simple convex problem. We also discussed a numerical 
implementation of a solver for this type of convex problem. 

There are several directions for further research. Some of 
them have been discussed in section III-CI and III-Di Other 
potential questions include: How can we investigate more 
properties of a given binary function by using fewer measure- 
ments? Is that possible to even characterize the motion and 
evolution of a shape by this method? 
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reconstructed signal 



Fig. 6. The ID signal is blurred with a Gaussian filter with cr = 5 generated 
by Matlab command fspecial. The standard deviation of Gaussian noise 
added is respectively 0.03, 0.05, 0.07, generated by Matlab command randn. 
The positions of the bars in the reconstruction are sometime not precisely 
equal to the original due to the present of the noise. The minor difference 
between the reconstruction and the true signal are highlighted by circles. The 
average computational time for reconstructions are respectively 0.05s, 0.04s, 
0.03s. 
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number of measurements 




0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.0 
noise level 



Fig. 7. The trade-off between number of measurements vs. the noise level 
obtained from empirical tests. The noise level is defined as the standard 
deviation of the gaussian noise. For each pair of fixed noise level and number 
of measurements, 1000 random tests are taken to get the average number of 
miss -identified Os and Is. In each test the signal length is 100 and the number 
of consecutive constant intervals is 10. The right figure shows the levelset 
curve of the left figure. Different curves correspond to different number of 
miss -identified Os and Is. 




### 



Fig. 8. The first image blurred with a Gaussian filter with a = 5, generated 
by the Matlab command f special. The amplitude of Gaussian noise is 
respectively 0.03, 0.05, 0.07, generated by the Matlab command randn. 
The average computational time for the reconstructions is respectively 0.81s, 
0.81s, 0.80s. 




Fig. 9. The second image blurred with a Gaussian filter with a = 5, 
generated by the Matlab command fspecial. The amplitude of Gaussian 
noise is respectively 0.1, 0.2, 0.3, generated by the Matlab command randn. 
The noisy image is on the left, the reconstruction is on the right. The average 
computational time for the reconstructions is respectively 0.04s, 0.03s, 0.03s. 




Fig. 10. The last image blurred with a Gaussian filter with a = 5, generated 
by the Matlab command fspecial. The amplitude of the Gaussian noise 
is respectively 0.03, 0.05, 0.07, generated by the Matlab command randn. 
The average computational time for the reconstructions is respectively 1.54s, 
1.20s, 0.98s. 



Appendix 



Theorem. 111.11 Assume uq is a binary solution of Auq = b. 
There exists no nonzero v G {Av = 0} such that 

{v{x) < 0, when uo{x) = 1 
v{x) > 0, when uo{x) = 

if and only if uq is the unique solution of (Pi), i.e. solving 
(Pi) recovers uq. 

Proof: If (Pi) has another solution other than uq, denoted 
as u\ then let v = u' — uq, we have Av = Auq — Au' = 0. 
Moreover, since u'{x) G [0,1], Vx, then when uq{x) = 0, 
v{x) — u'{x) — uo{x) > 0; when uo{x) = 1, v{x) = u\x) — 
uo{x) < 0, which contradicts the given condition on v. The 
other direction is similar. ■ 

Lemma. 111.41 Let T = [0, 1] where and 1 are identified, i.e. 
T = = {z : \z\ = 1}. Given 2n points on T who define 2n 
intervals on T, there exists a real trigonometric polynomial, 
whose spectrum is limited in [— n, n], vanishing only at those 
points and changes signs alternatively on those intervals. 

Proof: This conclusion is natural in the context of trigono- 
metric interpolation. We denote the 2n points by 



{Ck 



^}lZ,cS' = {z:\z\ = l}. 
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Let 



C 



2n 



k=i 



where C = 11^=1 <^/c^^^- Then u{z) can be written as 



u{z) = ^ 



akZ 



k= — n 



i.e. u is a trigonometric polynomial with spectrum limited in 
[-n,n]. Since {q} C C G for z G 5'^ we have 



2n 



2n 



/c=l 



therefore is real on S^. If we treat u{z) = u{e^^^^) as a 
function defined on [0, 1], it is easy to check that u vanishes 
on and only on {au} and -^u{e^^^^) does not vanish on {a/e}, 
and the conclusion follows. ■ 

Theorem. 111.51 If uq{x) is a 1-D binary signal that can be 
represented as in ^ with 2d consecutive intervals of ones and 
zeros, then by knowing the Fourier coefficients {au} for \k\ < 
d, we can recover uq through the convex problem (Pi). (Notice 
that uo{x) e R,Vx implies ak = cT-k^Mk, so essentially we 
only need to know {ak} for < k < d.) This result is optimal, 
i.e. precise reconstruction via solving (Pi) is impossible if 
knowing even less. 



Proof: By theorem III. 31 we only need to construct a 
discrete signal v = A^r] = J^~^{Sr]) that changes sign only at 
the endpoints of the intervals, where 5* is the sampling operator 
selecting {ak : | A:| < d}. Let the starting points of each interval 
be {si}, by lemma llL4l we can find a trigonometric polynomial 



Ed 



2-Kikt 



that changes sign only at { (s^ — ^) } . Let 



v{x) 



k=-d 



akC 



27rik-j 



XG{1,...,7V}, 



then V (or —v) satisfies the requirement. To show that this 
result is optimal, we only need to notice that a trigonometric 
polynomial with order less than d cannot have 2d zeros due 
to the fundamental algebraic theorem. ■ 



Theorem. 111.91 Assume u{x^y) is a 2D binary function 
with analytic jump curve and the average directional number 
of zero-crossings of u along the angle is denoted as 
Kq. If there exists a band-limited real function v{x^y) = 
Y.ij,k)en^jke'^'''^^''^^^^ defined on T^, where Q = {{j,k) : 
< d}, such that the jump set ofu{x, y) corresponds 
to the zero levelset of v{x^y), then Kq < 2d^\/0. 

Proof: In this proof we are discussing the problem on 
a torus T^, so every coordinates are automatically mod by 1 
without explicitly written to make the notations clear. 

Without loss of generality, we only prove the case when 
e (— 7r/4,7r/4]. The other case can be automatically proved 
by switching the x and y coordinates. 

Since the zero levelset of v{x^ y) are analytic curves, Ke{s) 
is a piece-wise constant function with respect to both s and 
9 with finitely many jumps, and it is easy to see that we can 



only prove the theorem for m {0 : tan6> G Q, 6> ^ 0, 1} 
because this set is dense in (— 7r/4, 7r/4]. 
For any 9 ^ {) s.t. tan 6> G Q, assume 



tan6> = p/q-, p < q e Z slyq coprime. 



Let 



Fs,e{t) = (t.s^ttdinO) modi, teR 

be a linear flow on T^, then Fs^e(f) is a periodic function 
with period q since Fs^e(f + ^) = (t ^ q^ s ^ (t ^ q) tan 0) = 
{t,s ^ttSinO) = Fs^e(t) modi. 

The main idea of the proof is based on the following 
observation: a whole period of Ps,6i can be split to q segments 
Lsm,o, m = 0, 1, 2, . . . , — 1. Therefore, counting the average 
intersection along the line segments Ls^e can be replaced by 
counting the intersection along F. The latter is easier because 
it can be reduced to counting the ID zero crossings of a 
trigonometric polynomial inside a period. 

Recall that 

Ls,e{t) = (t, 5 + nan i9) mod 1, s e [0,1], t G [0,1] 

We will first show that Fs^e{t)., t G [0,g] can be split to 
q segments Lg^^e, m = 0, 1, 2, . . . , g — 1. Indeed, when t G 
[m, m + 1] for m = 0, 1, 2, . . . , g — 1, it is easy to check 

Fs,e{t) Ls^^p/q^Q{t - m), 

so when t goes from to q, Fs^e(f) can be seen as the 
connected version of q line segments {i^s+mp/g,6> * ^ = 
0, . . . , — 1}. Since p and q are coprime, elementary number 
theory tells us that 

{mod {mp/q, l)}^=o,i,2,...,g-i = {mod {m/q, l)}^=o,i,2,...,g-i 

Similar to the notation of #Ls^e, we use #Ps,6> to denote the 
intersection of Fs^e^ t ^ [0, ^] with the jump set of i^o, then 

q-l q-1 
#Ps,6' = ^ #^s+mp/g,6' = ^ #^s+m/g,6' 



m=0 

Therefore, we have 



m=0 



q—1 ^mdr. 



Ko = COS J ^Lg^ods = cos j ^t^Lg^eds 




= cosO I ij^Fg^eds 
/q 



(15) 



That is to say, we replace #Lg^o in the definition of Kq by 
#Fg^Q. Now we start to evaluate #Fg^0. Since uo{x^y) is 
corresponding to the zero levelset of v{x^y), #Fg^Q is equal 
to the number of zero-crossings of v{x^y) along Fg^o, and we 
need to count the zero-crossings of v{Fg^o{t)) when t G [0, q]. 
Let v{t) = v{Fg^e{t)), recaU 



v{x,y)= ^ ajke 
ij,k)en 



2tt i{jx-\-ky) 
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we have 

ij,k)en 

so 

v{qt)= a^fee2^^^^e^^^^^'^+^^)^ (16) 

U,k)en 

Since Fs^e{t) is periodic with period q, so is v{t), thus v{qt) as 
a function of t is periodic with period 1. Since ^/p^-k^ < d 
in Vt, we have 



Then ([T6l) tells us that as a function of t can be 

expanded as a trigonometric polynomial with order no more 
than d^p^ + q^, therefore the zero-crossing of v{qt) is no 
more than 2d^p^ + q^, i.e. #Fs^e < 2d^p^ + q^. Plug it 
into ([T5]) , we have Kq < 2dco^6^p^ q^ /q = 2d. 



Theorem. IIII.ll If uq is the unique solution of (Pi), h = 
b -\- e is the corrupted measurement, u is the minimizer of the 
optimization problem 

min\\Au-h\\l s.t. < u < 1, 

u 

then when \\e\\ < h{A^Ouo) where h > is a small amount 
depending only on A and O^o (see details in the proof), 
Bin) = uq. 

Proof: We introduce the following lemma: 

Lemma: Let O be an orthant such that O D ker A = {0}. 
When \\e\\ < h{A,0) where h > is a small amount 
depending only on A and O, the solution w* of the linear 
programming problem 

min \\Aw — e\\ s.t. w e O 

w 

satisfies \w*\ < 1/2, Vi. 

At first we demonstrate that this lemma implies the con- 
clusion we need. Because uq is the unique solution of (Pi), 
theorem urn tells us that O^o H ker A = {0}. If u solves 

min \\Au-b\\l s.t. 0<u<l, 

u 

then because b = b -\- e = Auq + e, it is easy to see that 
= u — Uq solves 



min \\Aw 



s.t. w G ' 



From the above lemma, since ||e|| is small enough, we have 
\w^\ < 1/2, Vz. Therefore B{u) = uq and the conclusion 
follows. 

Now we go back to prove the lemma. 
Consider the linear programming problem 

min ||74i(; — e|| s.t. w e O 

w 

The duality problem is: 

max -/i^e s.t. A^/i G O, ||/i|| < 1 



The Karush-Kuhn-Tucker conditions of the optimal variables 
of the primal and duality problems are 

e O, A^/i* e O, ||/i*|| < 1 
Aw"" - e 



0,Vz 
if Aw* - e ^ 



(17) 
(18) 



(19) 



where the latter two are the complementary slackness con- 
ditions. From ([TS]) we see that {w*^A~^ ji*) = 0, so when 
A^* - e ^ 0, by ([19]) we get 

-e) = -e||(A^*,/i) 

= -e||(^*,A^/i*) =0. (20) 

Apparently (l2Ql) also holds when Aw* — e = 0, so it is satisfied 
anyway. Therefore, Pythagorean's theorem shows 

\\Aw*r = \\ef -\\Aw* -cf <\\4\ 

i.e. \\Aw*\\<M. 

By the alternative theorem [II.2I from O H ker A = {0} we 
know that there exists v = A^r] G int(O). Without loss of 
generality we assume ||7^|| = 1. Let be a positive number 
such that \vi\ = \{A'^r])i\ > 2/i,Vi, then h depends on A and 
O only. When ||e|| < h we have 

{w*,v) {w*,A'^r]) 
= (A^*,r^)<||A^*||||/.|| = 11^^*11 <||6||</.. 



Because w* e O and v e int(O), w*Vi > 0,Vi Therefore, 
w*Vi < ^iW*Vi < h,\/i. From \vi\ > 2h,\/i we have \w*\ < 
h/\vi\ < 1/2, Vz. That finishes the proof of the lemma. ■ 
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